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The quantization of the electromagnetic field in a three-dimensional inhomogeneous dielectric 
medium with losses is carried out in the framework of a damped-polariton model with an arbitrary 
spatial dependence of its parameters. The equations of motion for the canonical variables are solved 
explicitly by means of Laplace transformations for both positive and negative time. The dielectric 
susceptibility and the quantum noise-current density are identified in terms of the dynamical vari- 
ables and parameters of the model. The operators that diagonalize the Hamiltonian are found as 
linear combinations of the canonical variables, with coefficients depending on the electric suscepti- 
bility and the dielectric Green function. The complete time dependence of the electromagnetic field 
and of the dielectric polarization is determined. Our results provide a microscopic justification of the 
phenomenological quantization scheme for the electromagnetic field in inhomogeneous dielectrics. 
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I. INTRODUCTION 

The quantization of the electromagnetic field in a linear 
dielectric medium has been the subject of many investi- 
gations since the first treatment by Jauch and Watson 
In a homogeneous and non-dispersive medium the 
photon is associated with the transverse part of the field, 
which can neatly be distinguished from its longitudinal 
part. In contrast, in an inhomogeneous non-dispersive 
medium the transverse and the longitudinal degrees of 
freedom get coupled, which renders the generalization of 
the quantization scheme to that case less straightforward. 
However, the quantization can still be accomplished by 
employing a generalized transverse gauge, which depends 
on the dielectric constant 0, 0, 0, 0, 013 ■ 

For a dielectric medium with dispersion the quantiza- 
tion procedure described above runs into problems. In 
fact, since dispersion in a dielectric medium is inextrica- 
bly connected to extinction, one should account for the 
effect of losses in the quantization procedure. Huttner 
and Barnett HJsl were the first to use the Hopfield po- 
lariton model [l(| for the description of a dielectric with 
losses. To incorporate losses they coupled the dielectric 
polarization to a bath of oscillators, which causes a damp- 
ing of the polaritons. Subsequently, upon assuming the 
medium to be homogeneous, they were able to diago- 
nalize the Hamiltonian of this damped-polariton system, 
and to establish explicit formulas for the electromagnetic 
field and the dielectric polarization in terms of the diag- 
onalizing operators. As a direct application of the model 
they evaluated the change in the atomic decay due to a 
dielectric environment [ll| . Later on, their work has been 
reformulated and extended in various ways. In [l3 | and 
[l3| a simplified expression for the dielectric constant of 
the model was found. An alternative description of the 



model in terms of path integrals was given in [l^. In 
P*5j Laplace transformations were employed to simplify 
the diagonalization process considerably. Finally, in 
these transforms were used to formulate the equations 
governing the dynamics of the inhomogeneous model, al- 
though a complete diagonalization was not attempted in 
that case. 

As stated above, the treatment in [H is con- 
fined to the homogeneous damped-polariton model. For 
such systems, a systematic use of spatial Fourier trans- 
forms greatly helps in carrying out the diagonalization 
of the Hamiltonian. Somewhat later it was realized 
[H HI HI HHHmm that the quantization of the 
fields in inhomogeneous media could be achieved in posi- 
tion space by adding a noise term to the Maxwell equa- 
tions in a phenomenological way. In agreement with the 
fluctuation-dissipation theorem, one then postulates suit- 
able commutation relations for this noise term. Account- 
ing for the Kramers-Kronig relations of the dielectric me- 
dia is found to be essential in defining field operators with 
standard commutation relations. An alternative formula- 
tion p4j of the quantization procedure by means of aux- 
iliary fields has been shown to be completely equivalent 

m 

The phenomenological quantization scheme has been 
extended to magnetic and to anisotropic media, and to 
finite media with gain (see the reviews [2(| [27j and refer- 
ences therein). The scheme has been applied to atomic 
decay [2^, |2!J, |3(], 1^3, to energy transfer [32, and to 
resonant dipole interactions [33j. The phenomenologi- 
cal scheme has also been used to study the properties 
of electromagnetic field operators in systems with opti- 
cal cavities or beam splitters [341 1351 , and with dielectric 
slabs [3(j . For these optical components input-output re- 
lations have been derived. It was found that extinction 
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usually has adverse effects on nonclassical properties of 
light, such as squeezing, nonclassical correlations and en- 
tanglement |8^.l83|. 

Although the phenomenological quantization of the 
electromagnetic field in absorptive dielectrics has been 
very successful, its connection to the damped-polariton 
model has been established only for the special case of a 
homogeneous dielectric medium. In fact, one would like 
to see whether the crucial properties of the noise term, 
which are postulated in the phenomenological approach, 
could be derived from the polariton model in the gen- 
eral inho moge neous case, as has been noted on various 
occasions [22, EM EH • Since the model furnishes a pre- 
cise Hamiltonian description of the interaction between 
the damped dielectric and the electromagnetic field, one 
would have obtained in this way a microscopic justifica- 
tion of the phenomenological quantization scheme. To 
arrive at this goal, one has to express the noise term in 
the canonical variables of the model, so that its proper- 
ties can be determined. Once the noise term has been 
found, it can be used for the complete diagonalization 
of the Hamiltonian. In the following we shall show how 
this can be achieved by employing a Laplace-transform 
technique as in [l5llT^|. 

The paper is organized as follows. In Sec. |H] the 
model is defined and the equations of motion for the 
canonical variables are derived. The Laplace transforms 
of these are determined in Sec. IIHI As we shall see, it 
is essential to introduce forward and backward Laplace 
transforms for positive and negative time, respectively. 
Subsequently, in Sec. IIVI the Green functions of the 
inhomogeneous dielectric are employed to determine a 
space- and frequency-dependent source density, which 
is the analogue of the noise-current density in the 
phenomenological quantization scheme. The explicit 
form for this source density as a function of the canonical 
variables is derived in Sec. E] Once identified in terms 
of the canonical variables, some important properties 
of the source density can be derived, as presented in 
Sec. I VII In particular, we shall demonstrate that the 
damped-polariton model can be diagonalized in terms 
of the source density. The full time dependence of the 
electromagnetic field and of the dielectric polarization 
density can thus be established. The paper ends with a 
discussion and with some conclusions. 



II. EQUATIONS OF MOTION 

As a model for an absorptive dielectric interacting with 
the electromagnetic field we adopt the inhomogeneous 
damped-polariton system. In this model the polarization 
density is a continuous space-dependent variable. The 
damping is provided through the coupling to a bath of 
harmonic oscillators with a continuous range of eigenfre- 
quencies. The bath coupling constant depends both on 
the frequency and on the position. The electromagnetic 



field is coupled to the polarization according to the stan- 
dard minimal-coupling scheme. 

The Lagrangian density of the damped-polariton sys- 
tem is 0: 

C = \e Q E 2 - \p^B 2 + \pX 2 - \pulX 2 



duj u 2 Y 2 

dojVu X • Y u . (1) 



+ 1P I dLuY 2 -±p 
lo Jo 

-$V • (qX) - aA • X 



The electromagnetic field is described by the scalar po- 
tential <I>(r) and the vector potential A(r), with E = 
— V$ — A and B = V x A. We choose the Coulomb 
gauge in which the vector potential is purely transverse, 
so that V • A = 0. The dielectric degrees of freedom are 
described by a space-dependent harmonic variable X(r), 
with an associated eigenfrequency wo(r) and a density 
p(r). The polarization density is given by -aX, with 
a(r) a space-dependent proportionality constant. In the 
Coulomb gauge the scalar potential is a material vari- 
able that is given as the solution of the Poisson equation 
A$ = — £q V • (aX), with suitable boundary conditions 
at infinity. The electromagnetic potentials <&, A and the 
dielectric variable X are coupled in the usual way, with 
V ■ (oX) the bound charge density and -aX the bound 
current density. Hence, a(r) gives the strength of the 
coupling between the electromagnetic fields (or the po- 
tentials) and X(r). Damping is introduced in the model 
through a continuum of harmonic-oscillator bath vari- 
ables Y w (r), labeled by the frequency u>. The coupling 
to X(r) is determined by the bath coupling parameter 
u w (r). A schematic representation of the system and its 
independent parameters is: 



electromagnetic field A, II 



a(r) 



dielectric 



oscillator bath 



X,P p(r),c (r) 



Vu,(r) 



Y W ,Q W p(v 



Here II, P and Q w arc canonical momenta, which will 
be defined below. 

In writing we have used the same conventions as 
in 0. In particular, we have refrained from a rescaling 
of the physical variables. As is clear from (Q, the den- 
sity p could have been scaled away by redefining X and 
Y w . In this way we would have been left with the in- 
dependent (rescaled) coupling parameters a and (and 
the frequency ujq). Turning the argument the other way 
around, we could as well have chosen different density pa- 
rameters px and py multiplying the contributions with 
X and Y w . In that way a model would have been intro- 
duced that seems more general, but it is not. 
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Introducing the canonical momenta 

II = ™=£|)A, 

pX — qA, 



dA 

ax 

6C 

syZ 



= pY w - v u X, 



(2a) 
(2b) 
(2c) 



with 5 a functional derivative in the variable to, we find 
the Hamiltonian as 



H = 



dr 
1 



ln 2 + i(VxA) 2 

2e 2p V ' 



2p 
1 

a 



P 2 + \pulX' 
dtuQl + 



-A P 



o 



drdr 



— A z + - / dw^X-Q, 
, V • (qX) V • (a'X') 



8vre |r 



(3) 



We introduced the notation X' = X(r'), and likewise a' 
and V. The rcnormalized frequency a)o(r) is defined by 
o)q = ujq + J dcj-y 2 ,, where the integral is assumed 
to be finite at all positions r. 

To quantize the model we impose the usual commuta- 
tion relations 

[II(r),A(r')] = -ifi«5 T (r-r'), (4a) 
[P(r),X(r')] = -iH\5(r-r'), (4b) 
[Q w (r),Y^(r')] = -iH S(u - «/) I S(v - r'), (4c) 

while all other commutators of the canonical variables 
vanish. Here I is the three-dimensional unit tensor, while 
<5x(r) = I S(r) + VV(47rr)~ 1 is the transverse delta func- 
tion. 

With these 'primary' canonical commutation relations 
and the relations between the potentials and field opera- 
tors, one can derive the following 'secondary' commuta- 
tion relations 



[E(r),A(r')] = -6 T (r-r>), 

So 



(5a) 



[E 3 {v),B k {v')] 



SO 



e^V^(r-r'), (5b) 



where we used (12 all . Notice that these commutation re- 
lations are medium-independent, because in Eq. (0J we 
took the electromagnetic field to be canonically indepen- 
dent from the material variables X and P. 

In the Heisenberg picture the equations of motion for 
the canonical variables follow by evaluating the commu- 



tators with the Hamiltonian: 
A 

1 



In, 

So 



n 



x 



AA- 



/'o 

V 

P 



(P 



aA), 



P 



-pu>Q X [<xX]l 

£o 



l 

p j(, 



-pu 2 Y u 



(6a) 
(6b) 
(6c) 

AluVu Q w , (6d) 

(6e) 
(6f) 



In the second equation the subscript T denotes the trans- 
verse part of the vector, which is obtained by a convo- 
lution with the transverse delta function. Likewise, the 
subscript L in the fourth equation indicates the longitu- 
dinal part, obtained by convolving with the longitudinal 
delta function <5l(i") = — W(47rr) _1 . 
From lU)aJ|-||HcJ| one gets 



AA 



: A = po 



- (P + oA) 

P 



= p a [oX] T . (7) 



This equation is equivalent to Maxwell's equation 
-V x B + c~ 2 E = p Q aX, 



(8) 



since the Poisson equation for $ can be rewritten as 

[aX] L = -e V$. (9) 

The equations of motion (|6*c)l - (|6d|l for the variables of 
the dielectric medium yield a second-order differential 
equation for X: 

,oX + pujqX. = aA [aX] L / dwii^Qu, (10) 



so P Ja 

or alternatively, with the use of (|6*e|) : 



pX + /9cJqX = aA [oX]l 

s 



duj v u Y u 



(11) 



The change from ujq to luq should be noted here. The 
first two terms at the right-hand sides of l|l(J|l and 
are equal to -«E, as follows from JJjJ. 

Finally, the second-order equation of motion for the 
bath follows from (I6*cll-ll6fl) as 



(12) 



pYu, + pu> Y u = v u X. 



The second-order equations J7J), (|11|) and l|12|) determine 
the time evolution of the basic physical variables A, X 
and Y u , which represent the vector potential, the di- 
electric polarization density and the harmonic-oscillator 
displacement density of the bath. If the initial conditions 
of A, X and Y w and their first time derivatives are 
given at t = 0, these operators are known at any value 
of t, either positive or negative. 
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III. LAPLACE TRANSFORMS 

Our goal for the rest of this paper is twofold. First, 
we would like to determine the dynamics of the model 
and to find the complete time dependence of the canon- 
ical variables. Secondly, we want to show how the phe- 
nomenological theory emerges from our model. To that 
end, we must identify the elements of the phenomeno- 
logical theory (dielectric function, Green function, and 
noise-current density) in terms of the variables from the 
model. Then, after their identification we must show that 
these quantities have all the desired properties that are 
merely postulated in the phenomenological approach. 

The equations of motion of the previous section con- 
stitute a set of linear differential equations. In this sec- 
tion, we solve them in terms of their initial conditions, 
by introducing Laplace transforms. This technique has 
been used before in the analysis of the damped-polariton 
model |l5l llfij]. In particular, we will derive an equation 
for the Laplace transform of the electric field E, while 
taking care to eliminate the Laplace-transformed vari- 
ables describing the dielectric medium. To that end we 
shall first solve the Laplace-transformed equations for the 
bath variables. Subsequently, these will be used to ob- 
tain an identity that relates the Laplace transforms of 
the polarization density and the electric field. While es- 
tablishing that relationship, we shall identify the electric 
susceptibility in Laplace language. Finally, wave equa- 
tions for the Laplace transforms of the electric field will 
be deduced. All equations will be valid for an arbitrary 
spatial dependence of the variables of the model. 

For any time-dependent operator f2 the (forward) 
Laplace transform is defined as 



n(p) 



dt e - pt n(t) 



(13) 



Obviously, the Laplace transform contains all informa- 
tion on the time evolution of fl for positive t. In the 
following we wish to determine the time evolution of the 
relevant operators of our model for any time, either pos- 
itive or negative. Hence, we also introduce the backward 
Laplace transform: 



h(p) = / dte- pt Cl(-t). 



(14) 



Both transforms are defined for all p with Rep > 0. 

Carrying out the (forward) Laplace transformation of 
(|5e)l and lj6f|l. and eliminating Y u (p) we find 



QUp) = 



p 2 + OJ 2 

1 



Vuj X(p) 



. , [pQ,(0)-p W 2 Y w (0)] , (15) 

P z + LU Z 



with Q w (0) and Y u (0) the initial conditions at t = 0. Al- 
ternatively, we could have used l|12f) as a starting point. 
Upon performing its Laplace transformation, we may 



eliminate the initial condition for Y u with the help of 
l|pe|) . The initial condition of X is then found to drop 
out as well, so that l|15(l is recovered. 

From (|5cll - (j6d(l we find analogously, after elimination 
of P(p): 



1 



X(p) = - pA(p)--[aX(p)] L 



1 



dw» u Q u (p)+pX(0) + -P(0). 



Jo 

Insertion of i|15[l in the integral yields: 



P 



(16) 



P 



1 

7 



duj ■ 



- \pA(p) - -[oX(p)] L 
P I £o 



X(p) 
-pX(0) 



-P(0) 



1 

P Jo 



du> ■ 



P 



u 2 Y w (0)-^CL(0) 



P 



(17) 



The first few terms at the right-hand side are related to 
the Laplace transform of the electric field, which reads 

E(p) = - P A(P) + ~[«X(p)]l + A(0). (18) 

From (jT7Jl with 1(15)1 it follows that the Laplace- 
transformed polarization density -aX(p) is proportional 
to the Laplace-transformed electric field, apart from 
terms depending on the initial conditions. We identify 
the proportionality constant as the electric susceptibility 
in Laplace language: 



x(p) 



a 



1 



Jo/O p l 



,~,2 



P 2 f °° duu 2 v 2 /(p 2 +LJ 2 ) 



(19) 



All parameters at the right-hand side depend on position, 
so that the susceptibility is a space-dependent quantity. 
In this respect it is a generalization of the definitions 
m 

II I lH 0, which are valid for homogeneous di- 
electrics. The susceptibility is an analytic function of p 
for all p with Rep > 0. Indeed, one can prove that the 
denominator in l|19f) cannot vanish for any p in the right 
half-plane. As a consequence of its analyticity properties, 
the susceptibility satisfies the standard Kramers-Kronig 
relations which connect the real and imaginary parts of 

x(p)- 

After introduction of the susceptibility \{p) the linear 
relationship between X(p) and E(p) gets the form: 

X(p) = --x(p)E(p) 

a 



eo 



Xip) «A(0) + PP X(0) + P(0) 



du ■ 



P 



w 2 Y u (0)-^Q w (0) 
P 



(20) 



This identity, which relates the Laplace transforms of the 
polarization density and the electric field, will be cru- 
cial in eliminating the dielectric variables, as we shall 
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see presently. As expected, the relationship still depends 
on the initial values of all canonical variables (with the 
exception of 11(0)). 

Having succeeded in expressing the polarization den- 
sity in terms of the electric field, we want to establish 
a wave equation for the electric field in Laplace lan- 
guage. All Laplace-transformed equations, which we de- 
rived above, contain terms depending on the initial con- 
ditions. Hence, we expect that the wave equation will 
contain such terms as well. In fact, we shall show that 
owing to the presence of these terms, the wave equation 
will be inhomogeneous. 

To derive the wave equation in Laplace language, we 
have to use the remaining ficld-dcpcndcnt equations in 
the set ©. However, instead of l(6al) - (|6bl) we prefer to 
employ their corollary (JSJ, or, better still, the equation 
which follows by taking the time derivative of JBJ and 
using the induction law: 



V x (V x E) + c" 2 E = fi aX. 
After a (forward) Laplace transformation we get 



(21) 



V x [V x E(p)] + c" VE(p) - mo<V X(p) 
= c- 2 E(0) + c" 2 pE(0) - AioaX(O) - ^oapX(O). (22) 

As before, the introduction of Laplace transforms has led 
to terms depending on the initial conditions. We wish to 
express these in terms of the canonical variables at t = 0. 
At the right-hand side we use (|Sc|) for t = to rewrite 
X(0): 

X(0) = -A(0) + ip(0). (23) 
P P 

Upon inserting this relation in JHJ for ( = 0wc find 

2 

E(0) = c 2 V x [V x A(0)] + — A(0) + — P(0). (24) 



SoP 



Furthermore, we use l(6al) and @ to write the electric 
field at t = as 

E(0) = -— 11(0) + — [oX(0)] L . (25) 
So £o 

Substituting (23)- in the right-hand side of (|22J), 
and using the expression l|20|l for X(p) in the last term 
at the left-hand side, we finally arrive at the Laplace- 
transformed wave equation for the electric field in the 
form: 

V x [V x E(p)] + c"Ve(p)E(p) = J(p), (26) 

with e(p) — 1 + x(p) t ne (relative) electric permeability 
in Laplace language. 

The differential equation 1(26( 1 in Laplace language is 
an inhomogeneous wave equation. The source term at 



the right-hand side depends on the initial conditions at 
time t = 0: 

J(p) = --Vx[VxA(0)]-£„px(p)A(fl) + n(0) 
£ -^P 2 x(p)] X(0) - [oX(0)] L - -PX(P) P(0) 



1 



PXKP) / dw ^— 

a J Q p l + 



V 



WY u (0)-ilQ u (0) 
P 



.(27) 



In principle, the electric field in the inhomogeneous 
wave equation (12(i[| can be solved in terms of its source. 
Upon performing the inverse Laplace transformation, one 
then gets an expression for the electric field at any later 
time after t = 0. However, as noted above, we also need 
the electric field (and the other variables of the model) 
for all times prior to t = 0. To find that information 
we have to derive the corresponding wave equation for 
the backward Laplace transform of the electric field. We 
shall proceed along the same lines as above. 

The backward Laplace transforms of l|6*e|l and (|6T|l yield 
on a par with l|15|l : 



QUp) 



■ Vu; X(p) 



+ [ P Q u (0) + pw 2 Y u (0)] , (28) 

where the change of sign in the last term should be 
noted. Using this equation and following the same steps 
as above, we find the counterpart of (|20(l as 



X(p) 



so 



x(p) Hp) 



so 



Xip) -aA(0)+ppX(0)-P(0) 



da; 



p 2 + LU 2 



w 2 Y w (0) + ^Q w (0) 
P 



(29) 



The susceptibility is given by (|19|) . as before. The back- 
ward Laplace transform of the electric field is slightly 
different from (JTSJ: 



E(p) = P A(p) + - [oX(p)] L - A(0). (30) 

so 

Here, as in l|29|) . several terms have changed sign. Finally, 
turning to the field equation, we obtain from H21(l after a 
backward Laplace transformation: 

V x [V x E(p)] + c _2 p 2 E(p) - noap 2 X(p) 



; E(0) 



2 pE(0) + MoaX(0) -/ioapX(0).(31) 



Substitution of Q23 [l -I)25 [l yields as the analogue of 126(1 : 
V x [V x E(p)] + c- 2 p 2 e( P )E(p) = PopHp)- (32) 
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The source term in this inhomogeneous wave equation is: The adjoint equation of l|34[l reads 



J(p) V x [V x A(0)] - s px(p) A(0) - n(o) 

VoP 



1 - £ -^P 2 x(p) X(0) + [qX(0)]l - ^PX(P) P(0) 



— px(p) / dw 

" Jo 



p 2 + u 2 



W 2 Y^(0) + ^Q w (0) 
P 



(33) 



As expected, several terms have changed sign as com- 
pared to (|27l) . For future convenience we have chosen 
the sign in the right-hand side of ll-'il'l) to be the opposite 
of that in ffify. 

The main results in this section are the expressions 
(|2"T| and (|33|) for J(p) and J(p). These are the source 
terms of the Laplace-transformed wave equations Ij26|l 
and J221)- It should be noted that these source terms 
are not defined as the forward and backward Laplace 
transforms of some operator J(r, t), although their 
notation might suggest otherwise. As a consequence, 
their properties differ from other pairs of operators, like 
E(p) and E(p). In Sec. IVII we shall come back to this 
point. 



IV. GREEN FUNCTIONS AND SOLUTIONS OF 
WAVE EQUATIONS 



In the previous section we have seen that both the for- 
ward and the backward Laplace transform of the electric 
field satisfy a wave equation with a source term. To solve 
these equations we introduce tensorial Green functions in 
Laplace language. The Green function associated to the 
wave equations l|26|l and (|32l) is defined as the solution of 
the differential equation |22 l l23 l l2r | : 

v 2 

- V x [V x G(r,r',p)] - ^s(r,p) G(r,r',j>) 



G(r,r',p) x V' 



x V - ^■e(r',p)G(r,r',p) 

= I (5(r -!■'), (36) 



where the spatial derivatives operate to the left. 

In terms of the Green function, the solution of l|26|) for 
the Laplace transform of the electric field is 

E(r,p)= Mo p / dr'G(r,r»-J(r». (37) 



Likewise, the backward Laplace transform of the field 
follows from l|32l) as 



E(r,p) = -fM p J dr'G(r,r»-J(r». 



(38) 



IJ(r-r'), (34) 



The two integral representations l|37(l and (|38|l for the for- 
ward and the backward Laplace transform of the electric 
field contain all information that is needed to express the 
electric field at time t in terms of the initial conditions 
of the canonical variables. The latter show up explicitly 
when the expressions (|27|) and l-i-il) are substituted in the 
integrals. 

The time-dependent electric field is obtained from the 
integral representations by an inverse Laplace transfor- 
mation. From (|37|l we get the electric field for t > 0: 

/oo 
-oo 

x J dr'G(r,r / ,-ia; + 0) • J(r', -iw + 0). (39) 

Here we changed the integration variable from p in the 
right half-plane to — iu> + r), with a small but positive r\. 
We formally replace ij by 0, so that — iu + r\ becomes 
-id) + 0. 

The electric field for t < is obtained from the inverse 
Laplace transform of 138fl : 



here for clarity we reintroduced the spatial argument of 
the permeability. The function G(r,r',p) is the forward 
Laplace transform of the standard retarded Green func- 
tion of macroscopic electrodynamics. It is also equal to 
the backward Laplace transform of the advanced Green 
function. The Green function is analytic for all p in the 
half-plane with Rep > 0, as is the susceptibility x( r ;P) 
[23| . The defining equation l|34|l can be read as the state- 
ment that the Green function is the inverse of the opera- 
tor — [V x (Vx)] — c~ 2 p 2 e(r,p)(l-), which is a symmetric 
differential operator in the space of square- integrable vec- 
tor functions. Owing to the symmetry the Green function 
satisfies the reciprocity relation: 



[G(r,r»] = [G(v',r,p)] 



31 



(35) 



E(r,t) = / due-^u 



x J dr'G(r, r',iu + Q) ■ J(r', + 0).(40) 



The w-dependent integrand in l|39|l is analytic for all 
lj in the upper half-plane. Hence, the integral over u> 
vanishes for negative t, since the w-contour can be closed 
in the upper half-plane for t < 0. Likewise, the right- 
hand side of l|4UI) is zero for positive t. As a consequence, 
one may combine the two expressions into a single one, 
which is valid for all t: 



E(r,t) 



n OO 

/ due^'EW(r,w) +h.c, (41) 
Jo 
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with the positive-frequency Fourier component: 



Ef + I(r,w) = - 



IflOUl 

2tt 

dr' [G(r, r', -iw + 0) • J(r', -iui + 0) 

+G(r, r', ico + 0)- J(r', + 0)] . (42) 
Alternatively, one may write 



EW(i 



2tt 



E(r, -iw + 0) + E(r, iw + 0) , (43) 



as follows by going back to (f3"?jl and lf3*%|l . 

We have succeeded now in obtaining the full time de- 
pendence of the electric-field operator. As the representa- 
tion H41J) shows, the field is a linear superposition of con- 
tributions, each with its own time dependence, and with 
a weight that is determined by the positive-frequency 
Fourier component (|42|l . The latter is itself a linear com- 
bination of the canonical variables at t = 0, as follows 
from (J2T |) and ff3 |l . 

The representation (|41|) is valid for all t. Hence, con- 
sistency demands that the right-hand side of (|41|l should 
reduce to the electric field (|25|l in the limit t — > 0. In Ap- 
pendix A we show that this is indeed the case. It turns 
out that the proof depends on the validity of several sum 
rules for the Green function. The latter hold true as a 
consequence of the analyticity of the Green function and 
of its asymptotic behavior for large frequencies. 

Having determined the positive-frequency Fourier 
component of the electric-field operator, we may try and 
find the differential equation which it satisfies. Of course, 
one expects this equation to be of similar form as in i|26|) 
and i|32|) , with a frequency ui instead of the Laplace vari- 
able p. Let us introduce therefore the permeability in the 
frequency domain as e(r, a;) = e(r, — iui + 0), for real ui. 
Likewise, we will write x(r, — iui + 0) as x(r, ui). In terms 
of the permeability e(r, ui) we may define the differential 
operator — [V x (Vx)] + c~ 2 ui 2 e(r, ui)(\-). The positive- 
frequency Fourier component of the electric field satisfies 
a differential equation containing this operator: 



V x 



V x E(+)(r,u)l + ^e(r,w)E< + >(r,t 
= — ijiQU) J(r, ui). 



(44) 



As expected, the differential equation is inhomogeneous, 
like H26f) and (|32() . Apart from a trivial factor, the 
right-hand side contains a source term J(r, ui). Since 
the positive- frequency part E^ + ' (r, ui) of the electric field 
is known from l|42|l with l|27|) and (|33|) . the differential 
equation l|44|) may serve as the definition of the operator 
J(r, ui). Defining J(r, ui) in this way we try to make con- 
tact with the phcnomcnological quantization procedure, 
in which an equation of the same form as l|44l) plays an 
important role [H H E HJ EJ El . Although the 
equations have the same form in the two theories, their 
interpretation is rather different. In the phenomenologi- 
cal quantization procedure neither E^ + ^(r, ui) nor J(r, ui) 



are known at first. To obtain a well-defined theory one 
has to postulate several properties of J(r,w). In partic- 
ular, one postulates an identity for the commutator of 
J(r, ui) with its hermitian conjugate. Only after doing so 
does one arrive at a meaningful theory. In the present 
theory we are able to derive the algebraic properties of 
J(r, w), as we shall see later on. 

As in the phenomenological theory, the operator 
J(r, uj) may be interpreted as a (frequency-dependent) 
noise-current density. Its form still needs to be elabo- 
rated in more detail. By substituting (|42|l at the left- 
hand side of H44f> and using the definition (|34l) of the 
Green function to eliminate the differential operators, we 
get 



J(r,w) 



1 

2^ 



J(r, -iui + 0) + J(r, iui + 0) 



+ [X(r, -»w + 0) - x(r, iw + 0)] 

x / dr'G(r,r',!w + 0) • J(r',£w + 0). (45) 



It turns out that the operator J(r,w) has a more compli- 
cated structure than the positive-frequency Fourier com- 
ponent i|43|) of the electric field. Whereas the first two 
terms at the right-hand side have the expected form, an 
additional contribution shows up, which has arisen from 
the second term of (I42|) . In fact, the operator acting on 
the electric field in (|44|) is the inverse of G(r, r', —iui + 0), 
but not of G(r,r', iui + 0): the frequency arguments of 
the permeability do not match in the latter case. The 
additional term is proportional to the imaginary part of 
the susceptibility x( r , + 0). 

The operator J(r, w) will play an important role in 
the following. Before evaluating it explicitly in terms 
of the canonical variables, we will show that it can also 
be obtained in a different way. Let us consider, on a 
par with (|43|l . the positive- frequency Fourier component 
of X. As one might suppose that it is proportional to 
E( + ' , with a proportionality constant determined by the 
susceptibility, we will focus on the combination 



aXW(r,w 
a 



2tt 



£ox(r,w)E(+»(r, W ) 
X(r, -iui + 0) + X(r, iui + 0) 



r^xfr, -iui + 0) E(r, -iui + 0) + E(r, iui + 0) . (46) 

27T L J 

Comparing 122(1 and 126J1 and taking p = — iui + we infer 
that the forward Laplace transforms satisfy the identity: 

- aX.(-iui + 0) - e Q x{-iui + 0)E(— iui + 0) 



ui 
i 



-3{-iui + 0) 



-aX(0) +e o E(0) 



-aX(0) + e E(0)] , 



(47) 



where we suppressed the dependence on r for the mo- 
ment. Likewise, from l|31|) - (|32ll we derive for the back- 
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ward Laplace transforms: 



aX(iuj + 0) - sox{ioj + 0)E(iw + 0) 
= -3(ico + 0) + f-aX(O) + e o E(0) 

-i[-aX(0) + £ E(0)], 



(48) 



is due to a noise term in the polarization density, as 
has been noted for the homogeneous damped-polariton 
model 0, and in the context of the phenomenological 
quantization scheme [22I l2tj, [3(J . In the present model 
the noise-current density is a specific linear combination 
of the canonical variables, as we shall see below. 



where we note that the frequency argument in the suscep- 
tibility in front of E differs from that in the corresponding 
term in ljT7|) . Adding the right-hand sides of ilT7|) - (|4*g|l we 
see that all terms depending on the operators at t = 
drop out. The resulting equality can be used to evaluate 
the right-hand side of (|46() . if the susceptibility in front 
of E in l(4*5|) is changed to x(— iu> + 0) by hand. The 
correction term that is brought about in this way, can 
be rewritten by means of 138(1 . Finally, we arrive at the 
identity: 

-aXW(r,w) -eox(r,w)EW(r,w) = -J(r,w). (49) 

UJ 

Hence, the noise-current density J(r, uj) can also be 
found as the difference of the positive-frequency Fourier 
component of the polarization density — a(r)X(r, t) 
and £ox( r i Ci; ) times the positive-frequency Fourier 
component of the electric field E(r,t), apart from a 
trivial factor i/uj. Hence, the noise-current density 



V. EVALUATION OF THE NOISE-CURRENT 
DENSITY 

The expression (|45|l for the noise-current density 
J(r, to) is rather formal. It depends on the forward and 
the backward Laplace transforms J and J. In section iHll 
these have been given in terms of the canonical variables. 
By using (|27|l and l|33|) . we are able to express J(r,cj) in 
the canonical variables as well. 

We start by evaluating the first two contributions at 
the right-hand side of After substitution of (|27|) 

and l|33|) and adding the two contributions, several terms 
are found to drop out. The remaining terms are closely 
related. As before, we shall write x( r , —iuj + 0) as x( r 7 w ), 
and, correspondingly, X (r,iw + 0) as its complex conju- 
gate x*( r i ^-0- Furthermore, the imaginary part of x(r, w) 
will be written as Xi( r ; w )- Using these notations, we find 
from the first two terms of l|45|l : 



1 

2^ 



J(r, -iw + 0) + J(r, iuj + 0) = -f° w Xi (r, u) A(r, 0) + to 2 Xi (r, w ) X(r, 0) - ^ w Xi (r, J) P(r, 0) 

7T 



ira 



ira 



is 



do;'- 



w' 2 - (w + ioy 



- — wy*(r,w) / du' — 



lu' z - (w - iO) 2 



W ' 2 Y^(r,0) + -Q^(r,0) 
P 

c' 2 Y w ,(r,0) + -Q^(r,0) 
P 



(50) 



where all canonical variables are taken at t = 0. It should be noted that the two integral terms are not the hermitian 
conjugates of each other, since the signs of the terms with Q u > do not match. 

In the integral term of (|45|l we have to insert l|33(l for p — ito + 0. The term with the spatial derivatives of the 
vector potential does not drop out now, as it did in l|5U|) . It can be evaluated by a partial integration in r', which 
leads to an expression with derivatives acting on the Green function. Evaluating these by using (|36|l . we arrive at the 
following two terms: 



-^-uj 3 Xl (r,Lu) f dr , £*(r', W )G*(r,r', W )-A(r',0) + -wx l (r,^)A(r,0). 
7rcr J IT 



(51) 



Here we introduced the notation G(r,r',cj) = G(r, r', — iui + 0) (and hence G*(r, r',u>) = G(r,r',ioj + 0) as well), in 
analogy to the notations for e and \- The final term in (|51|l cancels the first term in l|5UI) . Part of the integral term 
in (|51|l (namely, with \* instead of £*) drops out as well, when the contribution from the second term in is taken 
into account. 

Collecting all terms, we arrive at the following result for the noise-current density: 



J(r,w) = J dr' {c A (r,r', W )-A(r , ,0) + Cn(r,r', W )-n(r',0) + c x (r,r', W ).X(r',0) + c P (r,r',^).P(r',0) 



+ / du)' Cyq(t,t' ,w,u>') 



/Y, ( (r',0) + 4Q u -(r',0) 
P 



(52) 
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The (tensorial) coefficients have the following form 



c Y Q(r,r',u;,u') 



TTCT 

c n (r,r',w) = — l — ui 2 Xi(r,cj) ^,(t,t',u>), 



7TC 



r,r ,lu) = w Xi( r ,w) \5(r - r ) -ui Xi( r , u) G T ,(r, r , u) 



cp(r,r',w) = 



-T7 W Xi(r,w)x (r , w) G (r,r,w), 

« Xi (r, W ) I <5(r - r') + W 3 Xi (r, w) x*(r', w) G*(r, r', w ), 

vr^, 7rc a 



+ 



7TQ 
£0 



Im 



X(r,w) 



I <J(r - r') 



u>' 2 - (w + iO) 2 

Xi(r,w)x*(r',w)G*(r, r',w). 



w 3 u'„ , 



Trc 2 a' uj' 2 - (uj - iO) 



(53a) 
(53b) 

(53c) 
(53d) 

(53e) 



In the first three formulas the complex conjugate of the 
Green function Gt' (r, r', w) appears. It is the transverse 
part of G(r, r',o;) with respect to r', which is defined 
by the convolution J dr" G(r, r", ui) ■ S^(r" — r'). In the 
last term of (I53cf) the symbol v' u , denotes the bath cou- 
pling parameter at the position r' (and the frequency 
ui'). Furthermore, a' and p' stand for a(r') and p(r'), 
respectively, as before. 

The coefficients (|53|) can be interpreted as commuta- 
tors. In fact, from (0J we infer 



~[J(r, W ),II(r',0)] = -c A (r,r», (54a) 



[J(r,w),A(r',0)] = Cn(r,r», (54b) 



[J(r,w),X(r , ,0)] = c P (r,r', 



(54c) 



VI. PROPERTIES OF THE NOISE-CURRENT 
DENSITY 

In this section we will determine a few of the properties 
of the noise-current density J(r, ui). In particular, we will 
focus on its commutation relations. 

We start by considering the commutator of J(r, ui) with 
the Hamiltonian To evaluate this commutator, we 
might use the expressions <|52|) - (|53|) . employ the canonical 
commutation relations (@J) and evaluate all contributions 
in a systematic way. Owing to the complexity of 15311 . 
this is a rather tedious task. A more convenient way to 
obtain the commutator is to use the expression l|45|) for 
J(r, ui) in terms of J(r, — iui + 0) and J(r, iui + 0). The 
commutators of the latter with the Hamiltonian can be 
found without difficulty. In fact, one gets by evaluat- 
ing the commutators of l|27|l and l|33(l with ©, or more 
straightforwardly, by employing the equations of motion 



and analogous relations for the other coefficients. 



It should be noted that the coefficients 153a|l - (|53djl are 
all proportional to the imaginary part Xi( r : w ) of the sus- 
ceptibility. Furthermore, the coefficient 153cfl is propor- 
tional to the bath coupling parameter v u . In the absence 
of absorption the dielectric is not coupled to a bath, so 
that v u vanishes. As (|19f) shows, the imaginary part of 
the susceptibility vanishes in that case as well, at least for 
all frequencies that are off-resonance. Hence, all coeffi- 
cients (I53fl are zero in this case, so that the noise-current 
density itself disappears. Clearly, the present formalism 
loses its meaning for a non-absorptive dielectric. 

As the noise-current density is fully known now in 
terms of the canonical variables, we can proceed and 
derive its properties. This will be the subject of the next 
section. 



and 



~[H,J(r,p)} = p3(r,p) 

— (v x [V x E(r,0)] + ^e(r,}>)E(r,0) } , (55) 



:[ff,J(r,p)] = -pJ(r,p) 



— ( V x [V x E(r, 0)] + V -e(v, p)E(r, 0) \ , (56) 
Hop I c z 



where the electric field is taken at time t — 0. Using these 
expressions, we find as the contribution from the first two 
terms in (|45|l to the commutator {i/K)[H, J(r, ui)\. 



iui 
2^ 



J(r, -iui + 0) + J(r, iui + 0) 
+ ^% ti (r,w)E(r,0). 



(57) 
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Furthermore, the last term in (|45(l contributes: 
-^ 3 Xl (r,w) / dr'G(r,r',iu + 0) ■ J(r',iu/ + 0) 



EO 



u;xi(r ) a;)E(r ) 0). (58) 



On adding the two contributions, we see that the terms 
depending on E(r, 0) cancel. The remaining terms are 
proportional to J(r, w), so that we arrive at the simple 
result 



j:[H, J(t,w) 



-ico3(r, lu) 



(59) 



To understand how this commutation property comes 
about, it is useful to give a somewhat more formal deriva- 
tion of the commutator. To that end we start by remark- 
ing that for an arbitrary operator £l(t) the commutator of 
the Hamiltonian with its Laplace transform Q(p) follows 
directly by Laplace-transforming the equation of motion 
in the Heisenberg picture: 



-[fr,fi( P )]=pfl( P )-fi(o). 



(60) 



Writing the analogous equation for the backward Laplace 
transform: 



-[ff,n(p)] = -pn(p) + n(o), 



(61) 



and adding the two equations after substitution of the 
appropriate arguments p, we get 



H, {Q(-iuj + 0) + n(ioj + 0)} 
-ioj{(l(-iuj + 0) + &(iui + 0)} 



(62) 



In particular, one gets for Sl(t) — E(r, t) by comparison 
with PI)): 



#,E (+) (r,w)l = -^E (+) (r, 



(63) 



Of course, this could not be otherwise: if it did not hold, 
the time dependence in (|41|l would be compromised. By 
invoking the definition l|44(l of the noise-current density 
in terms of the positive-frequency Fourier component of 
the electric field, it is immediately clear now that J(r, oj) 
must satisfy a commutation relation of the same form, 
which is indeed what we got in (|59|) . 

It should be noted that the commutator expression in 
(|55|l contains an additional term that differs from that 
in (|60|l . The reason for this discrepancy is that J(r,p) 
has not been defined as the Laplace transform of some 
operator J(r, t), as we noticed already in Sec. lIIII Similar 
remarks apply to J(r,p). 

Let us now turn our attention to the commutator of 
J(r, lo) with its hermitian conjugate at a different posi- 
tion and frequency. In view of the general form (|45|l it is 
convenient to start by calculating the three commutators 



involving J and J. These follow by substitution of i|27|) 
and and use of the canonical commutation relations 
(gj). The results are given in (|B1(I and l|B3(l of Appendix 
IbI As shown there, these commutators of J and J can 
be used to prove the commutation relation: 

J(r, W ),[J(r',a/)] 

= — w 2 X *(r, u) 5(u-u')l S(r-r'). (64) 

In an analogous fashion one may evaluate the commuta- 
tor of the noise-current density with its counterpart for 
different arguments. It is found that this commutator 
vanishes: 



[J(r,w),J(r', 



= 0. 



(65) 



As shown in appendix^] the commutators (|64|l and l|65|) 

appear as the results of calculations in which several 
terms cancel one another. In fact, all nonlocal terms 
involving transverse delta functions and Green functions 
drop out. The final answers show that J(r, uj) is a strictly 
local operator in its space variable: for all r' ^ r it com- 
mutes both with J(r',w) and with the hermitian conju- 
gate of the latter. Moreover, the right-hand sides of i|64|) - 
(|65|l show that the noise-current density is local in the 
frequency as well: the commutators vanish for u> ^ uj' . 
As a final comment we note that the commutator (|64f) 
is proportional to the imaginary part of the susceptibil- 
ity. Hence, it vanishes if there is no absorption. As we 
have seen above, the noise-current density itself vanishes 
in that case, so that l|64(l becomes a trivial identity. 

As demonstrated above, the commutator properties 
(|59|l . (|64|l and l|65|> follow from the dynamics of the 
damped-polariton model. These commutation relations 
are the same as the postulated relations of the noise- 
current density in the phenomenological quantization 
scheme [13, EE [T^ EH H3, H3- Evidently, the sta- 
tus of the commutation relations is rather different in 
both schemes. Our results provide a justification for the 
postulates in the phenomenological theory. 

The collection of operators J(r, uj) possesses another 
convenient property: together with their hermitian con- 
jugates they form a complete basis set for the canoni- 
cal variables of the model. It means that each of these 
variables can be written as a linear combination of the 
operators from the basis. To prove this statement, we 
may argue as follows. Let us tentatively write the vector 
potential as 

//■oc 
dr' J dujJ(r',oj)-f A (r',r,uj)+h.c, (66) 

with as yet unknown tensorial coefficients f^. Taking the 
commutator of both sides with [J(r", w')]' we find from 
(|54E|l . and P^): 



fA(r,r» 



r>), 



(67) 
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so that we get: 
A(r,0) 



dr' 



doj ■ 



1 



10 u 2 Xi{v'iu) 
xJ(r» •(£(]•', r,w) + h.c. . (68) 



To really establish the validity of this equality, which 



J 



we found by merely assuming the general form (|66|l . we 
insert H52[) in the right-hand side, which then becomes 
a linear combination of the canonical variables. Upon 
evaluating the resulting integrals with the techniques of 
Appendices and [51 we indeed find that only the term 
with the vector potential survives, and that the left-hand 
side is recovered. 



Two other examples of identities, which may be checked in an analogous way, are: 

- f dr' r duj - 

£o J Jo u; 2 Xi(r',uj) 

in f , , f°° , 1 

£o 



n(r,o) 

X(r,0) 



dr' 



du) 



■ J(r', uS) ■ c* a {y',y,uj) + h.c. 



J(r',w) • Cp(r',r,w) + h.c. 



(69a) 
(69b) 



Similar identities are found to be valid for the canonical variables P, Y w and Q w . Since all canonical variables can 
thus be expressed in terms of J and its hcrmitian conjugate, these operators must form a complete basis, as we set 
out to prove. 

The completeness of the set of operators J(r, u>) and [J(r, to)]', and their properties (|59|l . (|64|l and (|65|) imply that 
the noise-current density is proportional to the diagonalizing operator of the system Hamiltonian ©■ In fact, we may 
write: 



H 



IT 

£o 



(iv i duo 

o uxi{r,U) 



[JCr.wJl+.Jfrw) 



(70) 



We have checked that Q can be recovered from (|70|l . This is accomplished by substituting (|52[l with (|53fl in (|70|l and 
evaluating the resulting expression in terms of the canonical variables. A few details of this calculation are given in 
Appendix C. It should be noted that the two expressions for the Hamiltonian do not agree completely: they differ by 
a c-number, which corresponds to a zero-point energy. 

Now that we have succeeded in obtaining the diagonalizing operators of our model, we can determine the full 
time dependence of the vector potential, the electric field, the polarization density, or any of the dynamic variables 
that we have considered above. For example, the vector potential at time t follows from l|68|) by substituting the 
time-dependent noise-current density e~ luJt J(r', ui) in the integrand. Inserting the expression i|53b|) for Cn, and using 
we get 



A(r,i) 



-Mo 



dr' 



duj e- iut G T (r, r', w) ■ J(r', ui) + h.c. 



(71) 



In the phenomenological quantization scheme an integral representation of the same form shows up |22|. However, in 
that theory the noise-current density J(r, uj) is a formal operator. In the present model we have an explicit expression 
for J at our disposal. In fact, by substituting l|52|) we may evaluate the right-hand side of l|71|l in terms of the canonical 
variables at t = Q. The results are presented in Appendix[Dj As shown there, the vector potential gets a simple form 
when enough time has passed for transients to die out. In that long-time limit, it reduces to a linear combination of 
bath operators only: 



if i r°° 



VHx(r»G T (r,r» 



W 2 Y w (r',0) + -Q4r',0) 
P 



+ h.c. 



This expression has the same form as l|71f) . with the noise-current operator replaced by: 



J;(r,w) 



-7^v(oj)x(r,u) 



w 2 Y w (r,0) + — Q u (r,0) 



(72) 



(73) 



The combination between square brackets is proportional to the annihilation operator of the bath harmonic oscillators 
at the chosen position and frequency. Indeed, J; satisfies the same standard commutation relations l|6^l) - (|6*5|l as J. 
The vector potential thus depends on the bath annihilation and creation operators only, when all transients have died 
out. This result for the long-time limit is the generalization of a similar finding for the homogeneous damped-polariton 
model, which we discussed before [l5|. 
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The time-dependent electric field E(r, t) could be found in principle by separately evaluating its transverse part 
— n(r,t)/e from ((69a(l and its longitudinal part [aX(r,i)]j,/Eo from 169 bp . and adding the two contributions. How- 
ever, a simpler way to obtain E(r, t) is to insert the solution of 144JI into the general form 1(41(1 . In this way we get 

E(r, t) = -i/j,o J dr' J diu e~ lu,t wG(r, r', uj) • J(r', u) + h.c. . (74) 

In Appendix [D] it is shown that in the long-time limit the electric field is given by an expression of the same form 
as ((74(1 . with J replaced by J;. By using the commutation relations 1)64 (1 - 1)65(1 . which are valid for J/ as well, one may 
show that the commutator of the electric field and the vector potential in the long-time limit has the standard form 
(|5aj) . 

As a final example, we consider the time-dependent polarization density -aX(r,f). It follows from ((69b(l as: 

-aX(r, t) = -— [ dv' r dto L e-^ J(r', uj) ■ c* P (r', r, uj) + h.c. . (75) 
J Jo u 2 Xi{r',uj) 

Substituting ((53d(l we get: 

-aX(r, t) = — ~ dr' duj e^' uj x(r, uj) G(r, r', uj) ■ J(r', uj) + i / duo e~ iut - J(r, uj) + h.c. , (76) 
c J Jo Jo ^ 



where we used l(35() . This form for the time-dependent 
polarization density shows that it is the sum of a term 
involving the properties of the medium through the sus- 
ceptibility and a term which is determined by the noise- 
current density only. In fact, this is consistent with 
l(49(l . which was written in terms of the positive- frequency 
Fourier components. Indeed, the integrand in the first 
term at the right-hand side of ((76() is proportional to the 
positive-frequency part of the electric field, as we have 
seen in d}. In Appendix iDl the time-dependent polar- 
ization density is evaluated in terms of the canonical vari- 
ables at t = Q. Furthermore, it is shown there that the 
long-time limit of -aX(r, t) follows from ((76(l by replac- 
ing J by J;, as was found above for the vector potential 
and the electric field. 

The expressions l(7T|) . l(7^j) and ijTT))) give the com- 
plete time dependence of the vector potential, the electric 
field and the polarization density in the inhomogeneous 
damped-polariton model. For the special case of a homo- 
geneous medium, the expressions reduce to those given 
in 0. As we have seen, the implicit dependence on the 
canonical variables at t = can be made explicit by sub- 
stitution of the noise-current density in the form of ((52!) . 

Now that we have found in 171(1 and ((74(1 the explicit 
time dependence of the field operators A and E in 
terms of the noise-current density J, we can determine 
the commutators [E(r, t), A(r', t)] and [E(r, i), B(r', t)] 
for arbitrary t. Actually, for t = we have already 
determined these in JSJ and the fact that these commu- 
tators are medium-independent almost directly followed 
from the standard commutation relations (0J. Now that 
we have integrated out the dynamics of the material 
variables, the expressions ((71(1 and 1(74(1 clearly both 
depend on the medium through the Green function 
G(r, r',oj) and the noise-current density 3(r',u>). How- 



ever, since no approximations were made in order to 
obtain the time dependence of the field operators, their 
commutators should still be medium- independent, and 
equal to those at t = 0. Wi th t he use of ||21|)-(|2SJ| and 
the Green function sum rule l(A4(l one can can verify that 
the commutators [E(r, f), A(r', t)\ and [E(r, t), B(r', i)] 
indeed have the medium-independent values of |J5}. 
A medium-independent commutator [E(r, t), B(rVi)] 
was also found in the phenomenological scheme |23| . 
which was the principal argument in showing that the 
phenomenological scheme is consistent with (although 
not founded on) quantum electrodynamics. Finally, 
it may be remarked that in our theory the commu- 
tators [P(r,t),X(r',t)] and [Q w (r, t),Y u > (r',i)] are 
also medium- independent; these commutators have no 
counterparts in the phenomenological theory. 



VII. CONCLUSION AND DISCUSSION 

By solving the inhomogeneous damped-polariton 
model we have established a rigorous basis for the phe- 
nomenological quantization procedure, which has been 
used to describe quantum phenomena in linear lossy di- 
electrics with great success. Up to now such a firm basis 
was available for homogeneous dielectrics only, through 
the pioneering work of Huttner and Barnett || 01 • 

As a tool in our treatment we have used forward and 
backward Laplace transformations. With the help of 
these we solved the equations of motion for the canon- 
ical variables. The Laplace transforms of the electric 
field were shown to obey wave equations with source 
terms that could be expressed in terms of the canonical 
variables at time t = 0. Upon introducing the Green 
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function of these wave equations we were able to de- 
rive an expression for the positive- frequency Fourier com- 
ponent of the time-dependent electric field. The latter 
was found to satisfy a wave equation with a frequency- 
dependent source term that could be interpreted as a 
noise-current density for the inhomogeneous dampcd- 
polariton model. Explicit expressions for this noise- 
current density in terms of the canonical variables of 
the system have been derived. By establishing its alge- 
braic properties we could prove that it is proportional to 
the diagonalizing operator of the model. Once we have 
shown this, the timc-depcndcnce of all relevant opera- 
tors can be determined. As an illustration we gave the 
time-dependent expressions for the vector potential, the 
electric field and the polarization density. 

In order to show the internal consistency of our re- 
sults, we have derived and employed several frequency 
sum rules for the tensorial Green function and for the sus- 
ceptibility, namely $K ty -$K ty . {C6jl and {C8 )l -l|C9 jl . The 
outcomes solely depend on the high-frequency asymp- 
totic behavior of the Green function and the susceptibil- 
ity. In the present model this behavior is determined by 
the values of the (local) parameters a, p, loq and v^. The 
Green function sum rules are generalizations of 'veloc- 
ity sum rules' that have been derived for homogeneous 
dielectrics [9Hl5|. 

We have expressed all field operators in terms of the 
noise-current density operators J(r, lo), which also diago- 
nalize the Hamiltonian. These operators and their hermi- 
tian conjugates were proved to be local both in position 
and in frequency: any pair of them commute when taken 
at different positions and/or different frequencies. The 
locality in position is not self-evident a priori, as some 
of the canonical variables of the model, namely A and 
II, satisfy a commutation relation l|4a|l with a non-local 
transverse delta function. The positive-frequency Fourier 
component E^ + ) (r, lo) of the electric field is non-local in 
space as well: it does not commute with its hermitian 
conjugate at a position r' (and at the same frequency to). 
The locality of J(r,w) with respect to the frequency is 
connected to the validity of a generalized optical theorem 
(|B7ll for the Green function. When the independent fre- 
quency variables in this theorem are chosen to be equal, 
it reduces to the standard form of the optical theorem 

na 

The diagonalizing operators are not unique. For exam- 
ple, if one breaks up the noise-current density in terms 
of its canonical elements according to (|52|l , then for long 
times after the initial time t = one finds that the field 
operators are determined only by the initial bath opera- 
tors, since time-dependent coefficients of other canonical 
variables all decay exponentially fast. If only long times 
are considered, J(r,w) can be taken to be proportional 
to the initial annihilation operator of the bath harmonic 
oscillator at position r and frequency lo. We stressed this 
point in [l5l ] for homogeneous dielectrics. Other diag- 
onalizing operators can be constructed by transforming 
the noise-current density J(r, lo) with arbitrary unitary 



transformations U(r, r', lo), but these would not have the 
physical interpretation of noise-current density operators. 

Our solution provides detailed information on the dy- 
namical behavior of absorptive dielectrics. This informa- 
tion can be used to study dynamical processes like spon- 
taneous emission of guest atoms in inhomogeneous me- 
dia. For instance, transient effects in emission processes, 
which we studied before in homogeneous media |l5j . can 
now be investigated in the general inhomogeneous case. 
Local-field effects, which by their very nature are brought 
about by inhomogeneities in the medium, form another 
field of interest for which our solution may be helpful. 

To obtain our results we have employed a Laplace- 
transform technique that we used before Il5l. An 
alternative method, which was adopted in is based 
on a diagonalization procedure due to Fano ^(j- We 
have been able to carry out the diagonalization of the 
inhomogeneous damped-polariton model along those 
lines as well. Details of that work will be published 
elsewhere Eu. 
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APPENDIX A: SHORT-TIME LIMIT OF THE 
ELECTRIC FIELD 

In Sec. IL VI the electric-field operator E(r, t) has been 
found as a Fourier integral l|41|l . with the positive- 
frequency Fourier component E^ + ^(r, lo) given in l|42|l . It 
contains the source terms J and J, which have been given 
in H27|) and (|33|) as linear combinations of the canonical 
variables at t = 0. As a check, we shall verify that the 
Fourier-integral in H41[) reduces to E(r, 0) in the limit 
t — » 0. As we shall see, the proof will depend on the 
validity of a few sum rules for the Green function. 

By employing the identity [E(r, —ioj + O)]^ = E(r, iuo + 
0), and the corresponding identity for the backward 
Laplace transform, one may write the Fourier integral 
representing E(r,0) as 



h 1 ,u 



E(r, -iLo + 0) + E(r,-zw + 0) . (Al) 



Upon substituting dHZI> - <EHI> , 123 and lO, we find that 
all terms with A(r', 0), P(r', 0) and Y w /(r', 0) cancel. We 
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are left with the following expression: 



duj / dr' G(r, r', -iu + 0)- < Pquj n(r', 0) 



+/x w[a'X(r', 0)] T - + -^ 3 x(r', -*» + 0) X(r', 0) 



since the susceptibility x( r i ~^ UJ + 0) behaves like 
— [a 2 /(eop)]/(LO + iO) 2 for large w in the upper half-plane. 

Finally, the contribution of Q w is determined by the 
integral 



— 2-rX*(r ',-iw + O) 
x / cL/ ■ 



w' 2 - (w + i0) : 



Q^(r',0) . (A2) 



This result can be simplified by considering the integrals 
over ui. The terms with II and [aX]j contain the in- 
tegral did uj G(r, r', — + 0). Since the Green func- 
tion is analytic for u> in the upper half-plane, we may 
evaluate the integral by closing the contour in this half- 
plane. The Green function satisfies the differential equa- 
tion (|34|) . Now e(r,p) = 1 + x(r,p) tends to 1 for large p 
in the right half-plane, as follows from the expression l|19|) 
for x- Hence, the asymptotic form of the Green function 
for large ui in the upper half-plane has the same form as 
the free-space Green function, namely: 



dujuj 3 G(r,r', -iu+0)x(r', -iw + 0)- 



u' 2 -(oj + i0) 2 ' 
(A6) 

As the integrand is proportional to (uj + iO) 3 for large ui 
in the upper half-plane, a contour deformation leads to 
a vanishing result. 



From the above we conclude that the expression (|A2 
is equal to 



--n(r,0) + -[aX(r,0)] L , (A7) 



G(r,r',-iw + 0) 



(u> + iO) 2 



l<5(r-r'). 



(A3) 



As a consequence, closing the contour yields the identity: 



duo uj G(r, r', -iw + 0) = -inc 2 I 8(r — r'), (A4) 



which is in fact a sum rule for the Green function |2c 
Likewise, one proves 



duj uj 3 G(r, r', -icu + 0)x(r', -iui + 0) 



2 2 

vncra 



\S(r-r'), (A5) 



which is in agreement with (|25[) . 



APPENDIX B: COMMUTATORS OF THE 
NOISE-CURRENT DENSITY 



In this appendix we derive the commutators of the 
noise-current density with itself and with its hermitian 
conjugate. We start from (|45|l . in which J(r, u>) is given 
as a linear combination of the source terms J(r,p) and 
J(r,p). From the definition (|27|1 and the canonical com- 
mutation relations Q we obtain the commutator [J, «F] 
as 



[J(r,p), [J(r',p')] f ] = - V —^~ (W - I A) S(r - r') - fetf [ X (r,p) - x(r,p'*)] " 6(t - r'). (Bl) 

Mo PP p + p 



Here we used the auxiliary relation 



duj ■ 



2 2 



o (p 2 + uj 2 )(p' z + uj 2 ) 



a p 



£o p z — p 



x(r,p) x(r,p') 



(B2) 



where for brevity we did not write the position dependence of v^, p and a. The commutator [J(r,p), [J(r / ,p')]t] is 
equal to I |B1| I. apart from an overall minus sign. Finally, the commutation relation of J with «F is found to be 



J(r,p),[J(r',p 



' n'Mt 



— (W - I A) S(v - r') - * £o S [x(r,p) - x(r,j/')] I 8(v - r'). (B3) 

Mo PP P-P 



The last term at the right-hand side is not singular for p = p' . In fact, it is proportional to the derivative of the 
susceptibility in that case. It should be noted that the commutators (|B1|I and i|B3|> are local, as they vanish for r ^ r'. 
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This is not self-evident, since the definitions (|27|l and l|33|l contain nonlocal longitudinal terms. When evaluating the 
commutators, one finds that the nonlocal transverse delta function in the canonical commutator l|4a|l compensates for 
the nonlocal terms in J and J. 

Having derived the commutators for J and J, we can determine the commutator of J(r, uj) with [J(r', w')]t, by using 
Q45JI . In J we have to substitute p = —iuj + or p' — —id + 0, while in J the argument is p — iui + or p' — id + 0. 
Sorting out the various terms, we get three different types of contributions, namely those containing either no Green- 
function factor, or one or two such factors. In the first contribution the terms proportional to ( VV — I A) drop out. 
The terms proportional to I 5(r — r') cancel as well, unless uj and d are equal: 



1 

4^2 



J(r, -ico + 0) + J(r, iuj + 0), [J(r', -id + 0)] * + [j(r', id + 0) 



EqH 2 



d Xl (r,uj)5(uj-d)\S(r-r'). (B4) 



The contribution with a single Green-function factor is: 



e h 



iO 



Xi (r, u) Xi (r', d) d G(r, r', icj + 0) - d 2 G(r, r', 



0) 



(B5) 



Finally, the contribution with two Green-function factors is found as 
e h 



■ ud(u + uj 



+ 



eoh dd 3 



') Xi (r, oj) Xi (r', d) J dr" G(r, r", iW + 0) • { V" x [V x G(r", r', -id + 0)] } 
- x<(r, uj) Xl (r', d) { dr" [ X (r" , iu + 0) - x(r", -id + 0)] G(r, r", ioj + 0) • G(r", r', -id + 0). (B6) 



The second integral appearing here can be split into two parts, which may be rewritten with the use of (|34J) and 
After a partial integration we arrive at the identity 



dr" [ x (v", ioj + 0)- X (r", -id + 0)] G(r, r", ioj + Q) • G(r", r' 



0) 



doj' 2 



d G(r, r', ioj + 0) - d 2 G(r, r', -id + 0) 
+ {d - d 2 ) [ dr" G(r, r", ioj + 0) • { V" x [V" x G(r", r', -id + 0)] } 



(B7) 

For arbitrary uj and d this identity has the form of a generalized optical theorem for the Green function. By putting 
uj = d one recovers the optical theorem that has been discussed before [22j . When (|BT|t is used in i|B6|> , it turns out 
that all terms with spatial derivatives cancel, while the remaining terms are the opposite of (|B5|) . As a consequence, 
we are left with (|B4|) . so that we have proven the commutation relation l|64() . 

The commutator of J(r, uj) with J(r',ti/) can be evaluated in a similar way. As a preparation, one needs the 
commutators of J and J with their counterparts for different arguments. These operators satisfy the identities 



[J(r,p)]t=J(r )P *) , [J(r,p)]t=J(r,p*), 



(B8) 



as follows from inspection of (|27(l and (|33|l . Hence, the commutators of J and J can be written down immediately by 
using (|B1|I and l|B3|l . Subsequently, these commutation relations can be employed in evaluating the commutator of 
J(r, uj) with J(r', d). As before, we encounter terms with various numbers of Green-function factors. The contribution 
without a Green function is found to vanish: 



J(r, -iuj + 0) + J(r, iuj + 0), J(r', -id + 0) + J(r', id + 0) 



0. 



(B9) 



In an analogous way as above, the terms with one and with two Green-function factors in the commutator 
[J(r, uj), J(r',u/)] can be shown to cancel. This completes the proof of (JS3J). 



APPENDIX C: EVALUATION OF THE HAMILTONIAN 



In order to show how the Hamiltonian © is recovered from (|70|l . we give a few examples of the calculations that 
are involved. We shall demonstrate how the terms quadratic in II and in A in @ are obtained after substitution of 
(|52|l with (|53a|l - lj53b|l in l|7U|l . In this appendix all canonical variables are taken at t = 0. 
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We start with the term quadratic in II. After substitution of the appropriate expressions we get the following 
contribution to the Hamiltonian: 



1 



7reo c 



dr 



duj uj 3 Xi(r, oj) / dr' / dr" n(r') ■ G(r', r, -iu + 0) • G(r, r", iw + 0) ■ n(r") 



(CI) 



We could replace the transverse parts of the Green functions by the full Green functions, as II is purely transverse. 
Let us now rewrite Xi( r i w ) as — (*/2)[x( r : —i^> + 0) — x( r i + 0)] . Subsequently, we carry out the integral over r by 
means of the optical theorem, which follows from (|B7|I by taking u>' — uj. In doing so, the integral with the spatial 
derivatives in 1|B7J| drops out, whereas the contributions with a single Green function in (|B7|) remain. In this way we 
find from (|C1|) : 



dujuj dr' / dr"n(r')- [G(r',r",iw + 0) - G(r',r", -iw + 0)] • n(r"). 



(C2) 



The integral over the frequency can now be performed with the use of the sum rule (|A4|) . We finally obtain the simple 
result (2E0)- 1 / dr [II(r)] 2 , as in ©. 

As a second example, we consider the terms quadratic in A. From H7()J) with Ij52(l and (|53a|) we get a similar- 
expression as in (|C1|) . The main difference is a factor w 5 instead of w 3 in the integrand. After performing the integral 
over r as before, we arrive at 



te 



2?rc 2 J q 



duj uj 3 



dr' / dr" A(r') • [G(r', r", iuj + 0) - G(r', r", -iw + 0)] ■ A(r"). 



(C3) 



Owing to the presence of the factor w 3 we cannot use the sum rule (|A4j) . However, we may proceed by adding and 
subtracting the asymptotic form l)A3|l . In this way, the w-integral becomes : 



duj to 



G(r,r',-iw + 0) 



{uj + ioy 



\S(r-r') 



(C4) 



where we relabeled the position variables. Since the asymptotic form of the susceptibility for large uj is 
— [a 2 /(eqp)]/(uj + iO) 2 , as we have seen in appendix lAl one derives from the differential equation for the Green 
function the asymptotic form for large uj: 



G(r,r',-zw + 0) 



(uj + ioy 



IJ(r-r') 



(uj + ioy 



(VV- I A) S(r-r') 



Mo" 



I S(r - r') 



(C5) 



as a generalization of (|A3|> . Employing this asymptotic form in l|C4|l one finds by contour integration the sum rule: 



dujuj 



Substitution in (|C3|) yields 



G(r,r',-iuj + 0) 



(uj + ioy 



l<5(r-r') 



(VV - I A) 5{r-r') 



Mo" 



I5(r 



(C6) 



^ J dr A(r) ■ (VV - I A) ■ A(r) + J dr ^ [A(r)] 2 = ±- J dr [V x A(r)] 2 + j dr ^ [A(r)] 2 , (C7) 

which agrees with the contributions in Q- 

Similar techniques can be used to obtain the other terms in the Hamiltonian l(3|l. Several additional sum rules, 
with integrands containing products of the Green function and the susceptibility, are needed in establishing complete 
agreement. These are: 



dujuj 3 x(r, -iuj + 0)G(r,r',-iuj + Q)x(r',-iuj + 0) = 0, 
dui uj 5 x(r, -iw + 0) G(r, r', -iw + 0) x(r', -iuj + 0) = - 



l1TC a i xi i\ 
\5{r-r). 



Furthermore, one needs two sum rules for the susceptibility: 



duj uj x( r i ~i<- u + 0) 



mot 



duj uj 



x(r, -iuj + 0) + 



e a p (uj + i0) 2 



lira ujq 
£oP 



(C8a) 
(C8b) 

(C9a) 
(C9b) 



To prove these sum rules one uses a contour deformation and the asymptotic behavior of the integrands, as before. 
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APPENDIX D: TIME-DEPENDENT OPERATORS 



In this appendix we show how time-dependent operators can be expressed as linear combinations of the canonical 
variables. As examples we shall discuss the vector potential, the electric field and the polarization density. 

The time-dependent vector potential has been given in (|71|l . By substituting the expression (|52[l for the noise-current 
density, and employing the same methods as used in checking l|68|l we derive 

A(r,t) - J dr 'J ^e-^ t |c Aj4 (r,r',c l ;).A(r',0) + c An (r,r', W )-n(r',0) + c A x(r,r',c l ;).X(r',0) 



rOO 

■c A p(r,r',w) -P(r',0) + / duo' c AyQ (r, r', oj, oj') 

Jo 



W ' 2 Y^(r',0) + ^Q^(r',0) 
P 



h.c. 



The coefficients are 



c aa (y,y',oj 
CAn(r,r', oj 
CAx(r,r',w 
c ap {y,y',oj 



c AY q(y, t',oj,oj' 



1 



2 wlm [G T T'(r,r',w)] , 
Im [Gtt' ( r j r'; oj)] , 



7TC 

n 

7T 



Im [G T t' (r, r', oj)] oj 2 Im [G T (r, r', oj) x(r', 



-wlm [G T (r,r',w)x(r',w)], 



Im 



ojv'(oj') 



uo' z - (w + i0) : 



■ G T (r, r',oj)x(r',uj) 



(Dl) 

(D2a) 
(D2b) 
(D2c) 
(D2d) 
(D2e) 



The expression (|D1|I gives the vector potential A(r, t) for all t in terms of the canonical variables at t = 0. In 
particular, it may be used to determine the vector potential for large t, when all transients have died out. In order 
to derive this asymptotic form, one starts by noting that the term with the hermitian conjugate in IjDljl can be used 
to extend the w-integral over the whole real axis. For all positive t this integral may be evaluated by deforming the 
contour in the lower half-plane. If t gets large, the behavior of A(r,<) is dominated by the contributions from those 
singularities of the frequency-dependent integrand that are located close to the real w-axis in the lower half-plane. To 
find these singularities we consider the contributions from the various terms in (|D1|I one by one. The contributions 
involving A, II, X and P depend on the coefficients <|D2a|) - (ID2dll . The singularities in these coefficients are determined 
by those of (the analytical continuations of) x(r, oj) and G(y,y',oj). When the susceptibility has a finite imaginary 
part for real oj, continuity implies that the singularities of x an d G in the lower half-plane are located at a finite 
distance from the real w-axis. As a consequence, the contributions from the terms with A, II, X and P in IDlJl will 
die out exponentially fast for large t. On the other hand, some of the singularities of the coefficient (ID2e|) are really 
close to the real axis, as they are given by oj = ±oj' — iO. The contributions from these singularities will dominate 
the behavior of I|D1|I for large t. These contributions are readily evaluated by calculating the residues. One arrives at 
the result given in (|72|l of the main text. It shows that the long-time behavior of A(r,t) is governed by the specific 
combination J;(r, oj) of bath operators, given in JTSJ: 



A(r,t) 



-Mo 



dr' 



due- iut G T (r, r', oj) ■ 3i(r',oj) + h.c. 



(D3) 



As noted in the main text, the combination of bath operators occurring in J;(r, oj) is in fact proportional to the 
annihilation operator of the bath harmonic oscillator at r and with frequency oj. 

The timc-dcpcndcnt electric field E(r,i) can be evaluated in an analogous way. Substituting H52fl in H74f) we obtain 
an expression like (|D1(I . The coefficients Ce% (with i = A, II, X, P, YQ) follow from c A i in l|D2|) upon multiplying each 
coefficient by ioj and dropping the subscript T (but not T') of the Green functions. The analysis of the long-time 
behavior of the electric field is completely analogous to that of the vector potential. One finds an expression of the 
same form as (|73|) . with J replaced by J;: 



EM) 



-IflQ 



dr' 



doje- wt ojG(r,r',oj)-3i(r',oj)+h.c. 



(D4) 



Finally, we consider the time-dependent polarization density. From (|75|l and (|52|) with l|53|l we derive an expression 
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for X(r, t), which has the same form as (|D1|I . with the coefficients 



uj 2 Im [x(r, w) G T ' (r, r', ui)\ , 



Cxn(r,r» 
c M (r,r» 

c X p(r,r',w) 



1 



-wlm [x(r,cj)G T /(r, r',w)] , 



■ w Im [x(r, w) G T ' (r, r', u) ] 



W 3 Im [ X (r, w) G(r, r', w ) x (r', w)] , 



7ra 2 



Xi(r,w) l<5(r - r') 



«£0 



Im 



w' 2 - (w + iO) 2 
uj 2 «'(«') 



Im[x(r,w)G(r,r',w)x(r',w)] 
X(r,w) l(5(r-r') 



Im 



w' 2 - (a; + i0) : 



■X(r,w) G(r, r',w)x(r',w) 



(D5a) 
(D5b) 

(D5c) 
(D5d) 

(D5e) 



As before, the behavior in the long-time limit is dominated by the singularities close to the real frequency axis. These 
singularities arise from the coefficient in (|D5e(l . Evaluating their contributions we arrive at the following expression 
for X(r, t) in the long-time limit: 



X(r,f) 



dr' / duj e~ luJt u X {r, uj) G{r, r', uj) ■ J t {r', uj) 



1 

duje~ lujt - 3i{r,u) + h.c. , (D6) 

UJ 



with J; given by (|73|l . Upon multiplying with —a and comparing with we see that the long-time limit has the 
same effect on the polarization density as it has on the vector potential and the electric field: the dependence on the 
full noise-current density J(r, u>) is replaced by a dependence on J;(r,<x>). 
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